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■ A new Radiative Driven Implosion (RDI) model based on Smoothed Particle 

• Hydrodynamics (SPH) technique is developed and applied to investigate the mor- 

Q ■ phological evolutions of molecular clouds under the effect of ionising radiation. 

. This model self-consistently includes the self-gravity of the cloud in the hydro- 

^ . dynamical evolution, the UV radiation component in the radiation transferring 
equations, the relevant heating and cooling mechanisms in the energy evolution 

Q^. and a comprehensive chemical network. The simulation results reveal that under 

Q . the effect of ionising radiation, a molecular cloud may evolve through different 

. evolutionary sequences. Dependent on its initial gravitational state, the evolution 

^ci^ . of a molecular cloud does not necessarily follow a complete morphological evolu- 
tion sequence from type A^B^C, as described by previous RDI models. When 

^ . confronted with observations, the simulation results provide satisfactory physical 

0> , explanations for a series of puzzles derived from Bright-Rimmed Clouds(BRCs) 

\Q , observations. The consistency of the modelling results with observations shows 

^ ! that the self-gravity of a molecular cloud should not be neglected in any investi- 

^ ! gation on the dynamical evolution of molecular clouds when they are exposed to 

OO ! ionising radiation. 
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Introduction 



Bright-rimmed clouds (BRCs) found in and around HII regions are sites of ongoing star 
formation due to compression by ionisation/shock fronts, and provide an excellent labora- 
tory to study the influence of UV radiation from nearby massive stars on the evolutionary 
process of molecular clouds. Numerous molecular line, millimetre/sub-millimetre contin- 
uum and mid-IR surveys have revealed very detailed structures and physical properties of 



BRCs in various astro ph ysical environments (ISugitani et al 



Thompson et alll2004al Jbl: [Thompson fc Whitdl2004j : lUrquhart et al 



1991 



Sugitani fc Ogural Il994 
200d ). 



One of the most intriguing characteristics of the observed BRCs is the diversity of their 
morphologies. The observed BRCs appear in different morphologies even when they are in 
similar stellar environments. Sugitani et al (1991, 1994) classified 89 BRCs from a whole 
sky survey into three types depending on the curvature of the rim of BRCs. The three types 
are categorised with type A cloud with a rim displaying moderate curvature, type B cloud 
with a rim of a high degree of curvature which sometimes is also described as an elephant 
trunk morphology, and finally type C cloud with a tightly curved rim and a tail, which is 
also called comet ary globule. 



Lefloch fc Lazarefj (119941 ) presented a 2-dimensional numerical simulation of the effect of 



UV radiation on the dynamical evolution of a molecular cloud based on the Radiation Driven 
Implosion (RDI) model. Although their hydrodynamical modelling successfully created a 
complete sequence of the morphological evolution of a molecular cloud from type A to B 
and then to C, it seems difficult to provide physical explanations for the following questions 
derived from BRC observations: 

1) Why are there so many more type A BRCs? The first whole sky survey by Sugitani et 
al (1991,1994) revealed that 61% of the observed 89 BRCs have type A morphology. If all 
of the clouds followed the same evolutionary sequence revealed by Lofioch's modelling, from 
type A to type C morphology, there should be a balanced number of BRCs with different 
morphologies to be observed. 

2) Why do some of observed BRCs evolve to a quasi-stable type A rimmed morphology and 
show signs of star formation at the heads of their structures? If all of the clouds should evolve 
to type C morphology and finally were completely photoevaporated as described by Lelfioch's 
modelling, no star formation should be triggered in type A BRCs. The recent radio contin- 



uum and molecular line observations on BRCs in the catalo gue of SFO (ISugitani et allll991 



Sugitani fc Oguralll994l ) presented bv lUrguhart et all (120061) and many other observations in 



molecular lines and submilimeter lines (ILee et al 



2005 



Karr fc Martinll2003l : [Thompson et al 



2004bl ) provided strong evidences of ongoing star formations in the condensed head of type 
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A rimmed clouds where a pressure equilibrium at the boundary between the clouds and its 
environment has been reached. Therefore it seems that for some of the molecular clouds, 
further evolution to type B or C rimmed morphologies is not necessary, i.e., not every BRC's 
morphological evolution follows the sequence from type A to B then to C rimmed morphol- 
ogy, it may terminate their morphological evolution at any one of the three morphologies. 

3) What causes the spatial morphological distribution of BRCs in some of molecular cloud 
clust ers? The remnant molecular clouds in the Ori OB 1 associations (jOgura fc Sugitani 
19981 ) presented a spatial sequence from type A— i>B— i>C rimmed BRCs with their distance 
from the ionising star. It was once suggested that some r aolecular clouds which w ere farther 
away may directly evolve to type C rimmed morphology (lOgura fc Sugitanilll998l ). but there 
is no physical foundation for this speculation. 

For the molecular clumps very close to an ionising star, the ionised gas cannot expand 
freely from the surface of the clump, because their expansion is constrained by a material 
flow accompanying with the ionisation radiation flux. Lefioch's simulation showed that a 
zero radial veloc ity boundary condition, (the extreme case of flow effect) favours the type A 
BRC formation (ILefloch &: Lazarefllll994l ). This result may partly explain the predominance 
of type A BRCs at the border of the HII regions. However question remains for the type A 
BRC structures found at the outer location of the HII regions or typ e B/C BRCs found at 
the border of HII regions (ISugitani et allll99ll : ISugitani fc Oguralll994j ) and there are must be 
some other physical mechanisms which play dominant roles in the morphological evolution 
of molecular clouds. 

Therefore searching for satisfactory solutions to the above puzzles requires a more com- 
prehensive model which should include as many physical processes as possible so that reliable 
physical explanations on what observations revealed can be derived, hence to greatly improve 
our understanding on how the intensive UV radiation from nearby stars affects the dynamical 
evolution of the surrounding molecular clouds. 

The fact that IRAS sources and bipolar outflows in cloud globules are very often found 
at the edges of BRCs strongly suggests that the self-gravity of a molec ular cloud can play 
an ac tive role even in the early stages of evolution in molecular clouds (ILefloch &: Lazarefl 
19941 ). Nevertheless the self-gravity of a molecular cloud was neglected in Lefioch's model, 
because the initiative of the model was to investigate photo-evaporation effect of molecular 
clouds, for which self-gra vity does not play a dominant role. However, as pointed out by 
Lefioch fc Lazarefj (jl994j ) at the end of the discussions of their modelling results, inclusion of 
the self-gravitation of a molecular cloud is necessary in order to build a appropriate model for 
investigation on the dynamical evolution of globule clouds under the effect of UV radiation. 
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Recently, iKessel-Devnet fc BurkertI (120001 . l2003h have developed the first 3-D SPH RDI 



model for investigation of the effect of UV radiation on the dynamical evolution of molecular 
clouds. Their model includes both self-gravity and the hydrogen ionisation of the molecular 
cloud in the dynamical evolution equations and revealed further features of molecular cloud 
evolution under the effect of ionising radiation. However a self-consistent treatment for the 
energy evolution was not included which sets a barrier for a rigorous study on the evolution 
of the physical properties of a cloud, so that a direct confrontation with observations is 
difficult. 

Therefore it is our intention to build a comprehensive RDI model in order to obtain a 
consistent description about the effect of the UV radiation on the dynamical evolution of a 
molecular cloud and to reveal the physical origin for the observed characteristics of BRCs' 
morphology formation. 

In the following sections we first give an outline for our 3-D RDI model, then briefiy 
present an analytical solution which uses some observable physical properties of a condensed 
globule and the UV radiation flux to define a maximum mass a stable condensed globule could 
have to against the gravitational collapse. Next the dynamical evolution of the cloud will be 
detailed and the corresponding kinematics will be analysed based on numerical simulation 
results. Finally we will apply the derived knowledge to confront with BRCs' observations 
and the conclusions are drawn at the end of the paper. 



2. The model 

In this section, we present a brief description of the basic equations and various physical 
processes included in our model. The hydrodynamical equations are solved with Smoothed 
Particle Hydrodynamical(SPH) technique. For the deta iled description of the S PH theory 
and the corresponding technique, readers are directed to iNelson &: Langer 



2.1. The basic dynamical equations 

The continuity, momentum, and energy equations for a compressible fluid can be written 

as: 

| + pV-v = (1) 

^ = --VP-V<l> + Svisc (2) 
at p 
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-V- V 



7 



P P 
and the chemical rate equations take the general form: 

nKi 



dt 



(3) 



(4) 



where d/dt = d/dt + \ ■ V denotes the convective derivative, p is the density; v is the 
velocity; P is the pressure; Svisc represents the viscous forces; U is the internal energy per 
unit mass; and $ is the gravitational potential. 7 and A represent the non adiabatic heating 
and cooling functions respectively, which we will deal with later; The fractional abundance 
of main chemical species Xi in Equation (jlj) are for CO, CI, CII, HCO"*", 01, He"*", , OB.x, 
CB.X, M"*" and electrons; Ki is the associated chemical reaction rate and n is the total number 
density. 



2.2. Ionising radiation transfer 

The molecular cloud in our model is under the effects of both UV and FUV radiation 
fields. The UV radiation is from a nearby star and interacts with the molecular cloud through 
the top layer of the surface facing the star (front surface). The FUV radiation onto the front 
surfac e is from the nearby st ar and that onto the rear surface is from the background star 
light (INelson fc Langen Il999l ) . Then the FUV r adiation field around a mo lecular cloud can 
be approximated as a spherical radiation field (IGorti fc HoUenbachl |2002| ) . The the effect 
of the FUV field on the molecular cloud is mainly through the photoelectric emi s sion o f 
electrons from grains and is treated in the same way as that in iNelson &: Langerl (1l999l ). 
In the following part of this subsection we only need discuss how to deal with the ionising 
radiation (UV) field from the nearby star. 

Although Helium ionisation was included in the chemical network, we could safely ne- 
glect it when dealing with the ionisation radiation transfer, ionisation heating and cooling for 
simplicity, because of the much lower abundance of the Helium compared to hydrogen atoms 
( jPyson fc Williams! I1997I ). The implementation of ionising radiation from nearby stars into 
the above SPH code is based on solving the following ionisation rate and radiative transfer 
equations, 

due 
~dt 

— — = —an{l — x)J (6) 
dz 

where X = (jn{l — x)J is the ionisation rate with a being the effective ionisation cross section; 
TZ = nlas = x^ri^aB is the recombination rate with ub being the effective recombination 



i-n 



(5) 
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coeffi cient under the assumption of the 'on the spot' approximation (jPyson fc WiUiams 
19971 ). In its original definition, the recombination coefficient 



a 



E 



includes all of the individual recombination coefficients to the hydrogen atomic energy 
level i. In the 'on the spot' assumption, recombinations into the ground level (i = 1) do not 
lead to any net effect on the change in ionisation rate, since the photons released from this 
recombination process are able to re-ionise other hydrogen atoms on the spot. Therefore ai 
can be neglected and the resulting net recombination coefficient can be written as 



E 

i=2 



For a planner infall of ionising photons from a distant source on the boarder with a fiux Jq 
as the Lyman continuum photons per unit time and square area, the solution of Eq. ([6]) is 

J{z) = Joexp[-T{z)] 

where t{z) is the optical depth for ionising photons along the line of sight parallel to the infall 
direction of the photons, and z is the distance from the border of the integration volume 
along the line of sig ht. Neglecting the absorption e ffect of dust on ionising radiation, t{z) 
can be expressed as (jKessel-Deynet fc Burkertll2000l ) 



TiZ 



a dz' 



(7) 



R 



with hh = n{ l — x) and a being the mean o f over frequency, weighted by the spectrum 
of the source (jKessel-Deynet fc BurkertI l2000l ) . 



Eq.([6]) can be written in terms of the ionisation ratio x, i.e.. 



dx 
'dt 



1 dne 
n dt 



I 



n 



— X nas 



Solving Eq.® for the ionising ratio x includes three steps: i.e., calculate the optical path 
t{z); calculate the ion ising rate / and then fina ll y inte grate Eq.® for x. We employ the 
method introduced by iKessel-Deynet fc BurkertI (120001 ) in the above three-step work. We 
just give a ver y short description about i t here . The interested reader can find more detailed 
description in iKessel-Deynet fc BurkertI (120001 ). 



For any a simulated particle i with its number density n and its nearest neighbour 
list being given by the dynamical simulation at each time step, calculating t{z) is carried 
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through a) finding the evaluation points on the path toward source by smallest angel criteria 
and then b) summing up all the optical path elements on the path 

T{zi) = S™^o a nH,k 

where m is the total number of the evaluation points from the boundary of the volume to 
the position located one effective radius (oj = 0.85(Mj/pj)^/'^) before the point Zi which is 
derived from a); nH,k is the number density of the hydrogen atoms at the evaluation point 
k, Mi and pi are the mass and mass density of the particle i respectively. 

The ionisation rate for the particle i is then calculated by 

h = 7r^<^xp[-T{z,i)] [1 - exp{-2ai uha)] 

Finally, the first-order discretisation of Eq. (j8]) over a time interval At for particle i is given 
by 

^ = xl+At[^^- nl^' {x^^'an] (9) 
ni 

where the indexes j and j + 1 denote the values at the beginning and the end of the actual 
time-step 6t respectively. All of the values on the right-hand side are known from advancing 
the particles by the SPH formalism, except X/"*"^, which can be approximated by 

= n- — -%] (10) 



r - J+J- 
exp[—aai u 



2.3. The heating function 7 

The heating function in the thermal evolution is mainly provided by a) the photoelectric 
emission of electrons from grains 7pe caused by the incident FUV component (6.5 < hu < 13.6 
eV) of background star light, which can be treated as isotropic and has been discussed by 



Nelson fc Langeii (119991 ). so we don't repeat them here; b) the hydrogen ionisation heating 
'jion produced by the illumination of UV component {hu > 13.6eV) from a nearby massive 
star, which is highly directional and whose strength is dependent on the distance of the cloud 
from the star. Its formula can be derived through the following considerations. 

When an atom with ionisation energy Ei is ionised following the absorption of one 
photon of frequency u and energy E = hu, it releases one electron which carries the excess 
kinetic energy E^, — Ei which will be transfered to the gas through collisions with other gas 
particles. The resulting heating rate is described by the following equation. 

-fion = n{l - x)aJkT^ (11) 
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where /c is the Boltzmann constant and T* is expressed as (ICanto et allllQQSl ) 



xi + 2x0 + 2 



(12) 



where Tg// is the stellar temperature, and Xq = {Iiuq) j {kTej f) with z/q being the frequency 
of the Lyman hmit. 

Heating of the gas is also affected by cosmic ray heating, H2 formation heatirig and 
gas-dust thermal exchange, for which the same formulae as that in (INelson &: Langerlll999l ) 
are used in our model. 



2.4. The cooling function A 



The cooling function A is affected by recombination of the electrons with ions, the 
collisional excitation of Oil lines and by CO, CI, Cll and 01 line emission. 

When a free electron in the plasma is captured by a proton, a photon is emitted and an 
amount of energy Ei + mp v'^/2 is removed from the internal energy of the gas. The cooling 
rate due to this process is (IHummer fc Seatonlll963l ) 



Arec = PbIT- X kT 



(13) 



where /3b = x (1 + 0.158te 



Another important cooling resource comes from the collisional excitation of low-lying 
energy levels of Oil ions in spite of their low abundance, since Oil has energy levels with 
excitation potent ials of the order of kT. We use the following simplified formula derived by 
Raga et all (120021 ) to calculate the cooling rate due to the collisional excitation of Oil, 



A 



colli 



A, 



colli (1) 



+ A 



colli (2) 



with 



log 



10 [ 



log 



.A 



10 



rAcoiii(l), 

neUoii 

colli (2) 1 



7.9ti - 26.8 



1.9t 



10.5 



-20.5 



(14) 



(15) 
(16) 



neUon' \t2\ 

with ti = 1 — 2000 K/T and t2 = 1 — 5 x 10^ K/T and noii = xnoi- The above formula 
is valid under the limit of low electron density rie < 10"^ cm~'^, which is true for the BRCs' 
structure. 



The cooling a ffects Aime by CO, CI, Cll and 01 line emission are same as those in the 
Nelson fc Langerl Jl999[ ). 
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2.5. Boundary conditions and initial velocity field 



The molecule cloud is assumed to have a spherical symmetry and its mass is uniformly 
distributed through a sphere of radius R. It is embedded in su rrounding hot medium with a 
constant density (Hll region) (IKessel-Deynet &: Burkertll2003l ). The external pressure by the 
surrounding hot medium provides an inwardly directed pressure force whose characteristic 
magnitude is similar to that expected from a Hll region. The value of the assumed bound- 
ary pressu re Pe^i^t we adopted here is equivalent to an ionised gas with rig = lOOcm"^ and 
T = 10^K( lLefloch et alll2002l ). which resembles to Lefioch's zero radial velocity boundary 
condition. 

An outflow condition is imposed at a fixed radius equal to several times of the original 
radius of the cloud. Material that expand beyond this radius is simply removed from the 
calculation. By this time, these particles represent very low density material and so provide 
little UV extinction or other influence in the simulations. 

The initia l turb ulent velocity field was set up in the usual manner, as described by 



Mac Low et all (119981 ). with perturbations set up with a fiat power spectrum P{k), and with 
a minimum wavenumber kmin corresponding to a wavelength equal to the cloud diameter. 
Eight wave numbers were used to set up the velocity perturbations with the maximum, 
kmax = 8 kmin- The vclocity perturbations were drawn from a Gaussian random field de- 
termined by its power spectrum in Fourier space. For each three dimensional wavenumber 
k, we randomly select an amplitude from a Gaussian distribution centred on zero, and of 
height P{k), and selected a random phase uniformly distributed between and 27r. The field 
is then transformed back into real space to get the real velocity for each particle, and then 
multiplied by the amplitude needed to obtain the required total kinetic energy for the cloud. 
This procedure is repeated for each velocity component independently to get the full three 
dimensional velocity field. The amplitude of the turbulent motion was normalised such that 
the total kinetic energy of the clouds equals their initial gravitational potential energy. 

The initial geometry of the cloud and radiation field configurations are shown in Figure 
[U Since the distance of the cloud to the nearby illuminating star is assumed to be large 
compared to the radius of the cloud, the UV radiations from the star can then be treated 



as a plane-parallel radiation q r ito the front surface of the cloud (ILefloch &: Lazarefj Il994 



Kessel-Devnet fc BurkertI 



2000 



the cloud JGorti fc Hollenbachll2002l : iNelson fc Langei]ll999h 



Willams et al.ll200lf) while th e FUV radiation is isotropic to 



The radiation flux distribution at the boundary is 
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J(r = R, I 



Jo (Lyman) + Jo(FUV) 
Jo(FUV) 



0° < ^ < 90°) 
90" < ^ < 180°) 



where Jq (Lyman) is the value of the Lyman continu um flux at the surface of the cloud fo r 



which we take the typical value of 2 x 10^^ cm ^ s ^ f lMiao et alll2006l : IWillams et al.ll200ll ): 



Jo(FUV) is the value of the interstellar radiation at the surface of the cloud which has the 
value of 100 xGq) where Go is the rnean in terstellar medium radiation flux and traditionally 
known as the 'Habing flux' (IHabing Ill968l ); 6 is the azimuthal angle in a spherical coordinate. 



2.6. Classification of the initial condition of a cloud 



Lefloch &: Lazarefj (119941 ) classifled the initial physical properties of a cloud and UV ra- 
diation fle ld in a two dimensional p arameter space (A, F), which are deflned by the following 
equations (jLefloch fc Lazaret [l99J) , 



A 



uq 2r]aBnoR 
rjaBUiR 1 



(17) 
(18) 



where A relates the overpressure of the ionised gas to the undisturbed neutral gas and F 
is the ratio of recombination rate to the ionisation rate, = xn = Ug is the ionised gas 
density around the cloud, hq is the initial density of the hydrogen in the cloud, r] ^ 0.2 is a 
parameter to describe the effective thickness of the recombination layer around the molecular 
cloud and Cj = ^/VJT\ is the isothermal sound speed in the ionised gas with TZ being the gas 
constant. The 2-D parameter space can be divided into flve different regions according to 
the effect of the UV radiation on the evolution of a cloud. 

Region I is deflned by ■£ < 1.4 x 10^^, where the opacity of the cloud is very low 
( ^ 1 ) so that it is initially transparent to the UV radiation. Region II is deflned by 
A < ~ 10^^ — 10^'^ with Cn = yJTlTn being the isothermal sound speed in the neutral gas 
of temperature T„, in which the effect of the ionisation is too weak to produce any noticeable 
dynamical effect. Therefore clouds starting in region I and II are too trivial to be discussed. 

On the other end, when A is high, i.e., 2 < A < vTO, which deflnes the region III, 
the cloud is entirely photo-ionised in an ionisation flash by Lefloch's model in which the 
self-gravity of the cloud is not taken into account. Region IV is deflned by A < -^oiLyman) ^ 
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where the whole dynamical evolution of the cloud will be governed by the propagation of a 
D type ionisation front preceded by an isothermal shock in the cloud. Finally A > -^oiLynian) 
defines the region V, where an initial weak R ionisation front would gradually become a D 
type ionisation front propagating toward the interior of the cloud. 

For the clouds starting in regions III - V, although the detailed evolutionary features 
of the cloud can only be obtained through t he num e rical simulations described above, an 
existing simplified analytical result derived by iKahn I (119691 ) could help astronomers predict 
the final evolution of the globule under observation and also provides us a way to validate 
our numerical results. Therefore we briefly present it in the following subsection. 



3. A preview: A simplified analytical result 

In the simplified analytical result, only the effect of UV component is considered since it 
is much stronger than FUV radiation as given in section 2.5. The viscosity term in Equation 
(2) will be neglected because of its smaller magnitude than the other terms. 

After an initial ionising stage by the UV radiation flux, the stream of the ionised hy- 
drogen forms a warm halo around the globule which has a radius r^, at which the globule is 
bounded by an ionisation front while an isothermal shock propagates toward the interior of 
the globule. Assuming a steady ionised gas flow from the surface of the globule which can 
be approxi mated a s a sp herical isothermal sphere and a balanced ionisation-recombination 



in the halo, iKahn I (119691 ) derived a necessary condition for the globule to be stable against 
the gravitational collapse. Based on the analytical solutions of Equation (1-2) in the ionised 
gas halo and in the neutral globule, this critical condition is expressed by the current mass 
of the neutral globule 

M < M_ = 1.29( ^ , "if^ —t'^ (19) 

loTc^mifC^ Jo [Lyman) 

where G is the gravitational constant. 

The importance of this result is that the quantities such as T„(c„), Tj(cj) and Jo(Lyman) 
can be estimated by the observational data, then Equation ([19]) can be used to estimate the 
maximum mass for the globule to contain in order not to go collapse if the structures of both 
the globule and halo can be regarded as steady. On the other hand if the current mass M 
of the globule derived from the observational data is higher than Mmax, the globule under 
the investigation would finally collapse to form new star at the end of the evolution, vice 
verse, the globule may stay stable until it is completely evaporated by the UV radiation. 
The advantage of the Equation (|T9l) is that the criteria for globule's stability is in terms of 
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the current mass of the globule instead of the initial mass of the cloud, the latter is hardly 
known by observations. 

In the following, we will discuss the detailed evolutionary features of clouds, based on 
the numerical simulations. In this paper we will mainly study the morphological evolution 
of BRCs and try to find solutions for the questions raised in the beginning of the paper. A 
systematic investigation on the ionisation radiation triggered star formation in BRCs will be 
presented in our next paper. 

4. Results and discussions 

The most important difference between our newly developed RDI model and Lefloch's 
model is the inclusion of the self-gravity of a molecular cloud, therefore we will firstly in- 
vestigate the role of the self-gravity of a molecular cloud in its dynamical evolution when 
it is under the effect of UV radiations, so that we can explore the possibility to relax the 
extreme zero velocity boundary condition in order to derive the predominance of type A 
BRCs, and discuss the possible modification to the classification of Lefioch's 2-D parameter 
space caused by the inclusion of the self-gravity of the molecular cloud. 

We then will explore the effects of initial thermal state of a cloud, the strength of UV 
radiation flux, and also initial turbulent state (non-thermal motion) on the morphological 
evolution of molecular clouds. 

The molecular cloud is presented by 20,000 SPH particles. In the following analysis to 
the numerical results, we will concentrate on the effect of the UV radiation falling onto the 
surface of the front hemisphere, since the heating effect of photoelectric emission of electrons 
from dust grains by the isotropic FUV radiation is much weaker than that of hydrogen 
ionisation by Lyman continuum radiation which only affects those gas particles at the front 
surface layer. 

4.1. The effect of initial self-gravity of a moleculcir cloud 

The gravitational potential energy of a spherically uniform molecular cloud of mass M 
and radius R is Vt = — In order to explore the role of the initial self-gravity of a 

molecular cloud in its morphological evolution, we perform two different sets of simulations, 
one set of the three clouds having a same mass of 35 but different radii, the other having 
a same radius of 0.5 pc but different masses. Therefore our simulations cover different cases 
in which the initial self-gravity of a molecular cloud changes. 
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The initial temperature of the molecular clouds is all set to 60 K. These chosen initial 
conditions guarantee that the simulated molecular cloud will not collapse by self-gravity if 
there is no external radiations. 



4- 1.1. Clouds with a fixed mass of 35 Mq 

Due to the change of the radii of clouds, the initial ionisation/recombination states 
change as well. The first six columns in Table 1. list the parameters of the first set clouds 
called A, B, C and their initial locations in the ( A, F ) parameter space respectively. The 
features of cloud D in the list will be discussed later. The ratio of the initial gravitational 
energy of the three cloud A, B, and C is Qa '■ '■ = 1 '■ 0.38 : 0.26 and the ratio 
of gravitational forces of the clouds A, B and C onto a mass element at their surfaces are 
correspondingly Fa '■ Fb '■ Fc = 1 : 0.14 : 0.07, which defines significantly different initial 
self-gravity in the three molecular clouds. 

According to Lefioch's model, cloud A would evolve into a cometary globule after passing 
a morphological sequence from type A to type B then to type C, and finally it will be 
completely evaporated in a few Mys, while cloud B and C would be evaporated in a so called 
'ionisation flash'. It is our interest to see how these clouds would evolve differently if the 
self-gravity of a cloud is taken into account. 

I: CLOUD A 

The images in Figure [2] reveal the morphological evolution of cloud A by number density 
of hydrogen nuclei and the line plots in Figure [3] describe how the UV radiation induced shock 
front leads the ionisation front propagating into the interior of the cloud and finally results 
in the collapse of the cloud. Figure H] is the corresponding velocity field evolution diagram, 
which gives the kinematic explanation on the formation of the morphologies of the cloud at 
different evolutionary stages. 

The image in the top-left panel of Figure [2] and the solid line in the top-left panel of 
Figure [3] show the number density distribution of cloud A at an early stage (t = 700 years). 
It is seen that the density profile of the whole cloud has not yet changed much from the 
initial uniform distribution. The top-left panel in Figure H] shows a rather random velocity 
distribution which is very similar to the initial random velocity distribution of the cloud. 

However the dashed line in the top-left panel of Figure [3] tells that an ionisation front 
has been built up at the front surface of the upper hemisphere of the cloud after the intensive 
Lyman continuum flux falls onto the front surface. The gas in the top layer is fully ionised 
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[x ~ 1) and the resulted ionisation heating rises its temperature to lO^K. When t = 0.2 x 10^ 
years, the ionised gas at the front surface of the cloud evaporate radially away from the front 
surface as shown by the velocity profile in the middle panel of Figure H] and by the furry gas 
above the front surface in the top- middle panel of Figure [2J The falling tail beyond z = 0.4pc 
in the solid line of top-middle panel in Figure [3] shows that the density of the evaporated 
gas is about 100 cm~^. 

On the other hand, the high temperature in the front top layer produces a high pressure, 
i.e., an isothermal shock driven into the cloud. The neutral gas ahead the shock front is then 
compressed. The contour lines around the front surface layer of the cloud in the top-middle 
panel in Figure [2] and the small peak (if we define it as the shock front) around z = 0.4pc 
in the solid line of the top-middle panel in Figure [3] show that a density gradient starts to 
built up as a consequence of the shock compression. 

The UV radiation induced sh ock velocity Vc{6) at the front surface is approximately 
given by lLefloch fc Lazareffl (jl994j ). 



vM = vMicos{e)y/' (20) 

where 6 is the angle of a surface point from z axis and O<0<7r/2as shown in Figure 
[T] V^(0) is the shock velocity at the the point r = R,9 = and V^(0) ~ [Jo (Lyman) /n]^/^ 



(lBertoldilll989l ). The direction of Vc is along the normal of the point (R,^^). The resulted 
velocity distribution Vc{9) in the shocked gas is obviously seen in the top-middle panel of 
Figure HI 

At time t = 2.0 x 10^ years, the shock front leads the ionisation front moving toward 
the rear hemisphere as shown by the solid and dashed lines in the top-right panel of Figure 
[31 The ionisation ratio x has decreased a lot in the shocked part of the cloud, because the 
recombinations of electrons with their parent ions in the surrounding envelop have consumed 
most of Lyman photons. On the other hand the shocked gas at the front surface around 
^ ~ have the highest Vc so it moves faster toward the rear hemisphere than that in the 
peripheral part of the compressed layer. Consequently, the faster progression of the central 
parts of the compressed layer at the front surface causes a deviation of the morphology 
of the front surface from a hemisphere, i.e, the front surface of the hemisphere becomes 
seriously squashed so that a type A rim forms, which can been seen from the top-right and 
the subsequent panels of Figure [2l 

At the same time, the strongest compression of the gas located at ~ forms a small 
but highly condensed core at the head of the front hemisphere, shown as in the top-right 
panels of Figures [2] and [3l which creates a new gravitational centre G to its surrounding gas 
particles so that the latter gain an acceleration ac hence a velocity component Vq as shown 
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in Figure [51 The surrounding gas then gradually collapse toward the core centred at G, as 
shown in the panels from the top-right to the bottom-right panels in Figures [3] and |H 

Due to the distortion of the front surface, the velocity of the shocked gas in the outermost 
parts of the front hemisphere ( ~ 7r/2 ) is significantly changed so that its trajectory 
intersects the rear edge of the original cloud to form the 'ear-like' structures as shown by 
the two convex parts at ~ 7r/2 in the top-right and bottom-left panels of Figure O 

As illustrated in Figure \5\, the total velocity Vr of the gas particle at the sample point 
in the 'ear-like' structure directs to the C point on the z axis, instead of C point, which 
would be the converging point if the self-gravity of the globule was not taken into account. 
Hence the 'ear-like' structures converges to z axis to form a tail structure which overshoots 
the rear surface of the cloud and is well described by the solid line beyond z < — 0.5pc from 
the top-right to bottom-right panels of Figures [3] and HI 

At t = 0.351 My, the type A BRC structure comes to a thermally quasi-static state (with 
a pressure P/fc ~ 10^ Kcm~^ close to the boundary between the neutral gas and the ionised 
envelop) after an initial drastic and transient evolutionary period. As shown in Figure El 



the w hole structure has shown the typical stratification of an observed BRC (ILefioch et al 



2OO2I ) from the centre to outside: 



i) An isothermal neutral core of a radius of 0.15 pc (—0.3 < 2; < pc) with a central 
density of a few 10^ cm~^ and an average temperature of T„ = 28 K 

ii) A very thin layer of PDR (Photon-Dominated Region: < 2; < 0.05 pc) with an 
average density 10'^ cm~^, an temperature ~ 1000 K and an electron density Ue < 10 cm~^ 
with the ionisation fraction 10~^ < x < 10~^. 

The mass included in the isothermal core including the PDR (with density from lO'^ to 
a few 10^ cm~^ as observed) is 15 Mq. 

iii) A photoionised envelope (Bright rim: 0.05 < z < O.llpc), the density drops sharply 
from a few 10^ to a few 10 cm~^ and the ionisation ratio 10^^ < a; < 1, i.e., rie distribution 
between 10 and lO^cm"'^. The average temperature Tj is about 10^ K. 

iv) The evaporated ionised tenuous gas into the HII region {z > O.llpc) with an average 
density a few 10 cm~^ and x = 1, but a very high average temperature of 5 x 10^ K. 

Now the most concerned question on its further evolution is that whether the whole 
globule is stable against the gravitational collapse or not. Before revealing the numerical 
solution to this question, we can employ Kahn's analytical result to predict this globule's 
final evolution and then use the result to validate the numerical solution. Substituting all of 
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the values into Equation f[T9|) . we obtain the maximum mass for the isothermal core to be 
stable is M^ax = 3 Mq. Therefore this type A BRC should finally collapse since the current 
mass of the core is 15 Mq, much higher than M^ax- 

The numerical results after t = 0.35 My as shown by the bottom middle/right panels 
in Figures [21 [3] and H] reveal that the globule quickly collapse toward its G centre under 
the effect of the enhanced self-gravity. The central density of the core reached 2xl0^cm~^ 
when t = 0.391My and 5x10^° cm"'^ when t = 0.394My. Therefore we could conclude that 
the radiation triggered star formation in cloud A is expected to occur because the central 



density has reached such a high value (INelson fc Langerlll999l ). which is consistent with the 
analytical result. 

II: CLOUD B 

The morphological evolution of cloud B is shown in Figure [6l while Figure [7] describes 
how the shock front leads the ionisation front propagating into the cloud and Figure [8] 
displays the corresponding evolutionary sequence of its velocity field. 

During the first 0.4 My, the main characteristics of the morphological evolution of cloud 
B are similar to that of the cloud A. The gas particles at the top layer of the front hemisphere 
are continuously ionised, evaporated and form a hot photo-ionised envelope around the type 
A rim at the front surface of the globule. On the other hand, the radiation induced shock 
leading an ionisation front propagates into the cloud and compresses the neutral gas ahead 
to form a condensed core. 

However the morphological evolution of cloud B starts to show obviously different feature 
at t = 0.63My. It is seen from the bottom-left panel of Figures that the front surface 
becomes a type B rim with a dense core behind the rim. Also a more extensive tail than 
that in cloud A is formed by the gas particles from the 'ear-like' structure converging to z 
axis as shown by the solid lines beyond z < — 1 pc in the bottom panels of Figure [3 We will 
discuss the kinematics which make the difference in the morphological evolutions between 
cloud A and cloud B after the presentation of the evolutionary features of cloud C. 

Around the front hemisphere of the type B BRC, a similar thermally quasi-stable strat- 
ification to that of cloud A at t = 0.351 My is formed with an isothermal condensed core 
of a radius of 0.2 pc at the head of the globule as shown in the bottom-left panel of Figure 
[71 The core has a central density slightly lower than 10^ cm~'^, and an average temperature 
Tn = 31 K, which is surrounded by a thin layer of photo-ionised envelope of a temperature 
Tj = 10^ K. The core including the PDR contains a mass about 13.5 Mq which is higher 
than the Mmax =3.6 M0 given by Equation f[T9|) . Therefore the core should be expected to 
collapse. The simulation results presented in the bottom middle/right panels of Figure[6land 



-17- 



Figure [7] show that the core continuously gets denser, and the central density has reached 
a few xlO^ cm~^ at t = 0.79 My and a few xlO^° cm~^ at t = 0.81 My, the sign of the 
triggered star formation which is consistent with the analytical prediction. 

Ill: CLOUD C 

The evolutionary sequences of density, ionisation ratio and the velocity field of cloud C 
are shown in Figures [TUl and [TT] respectively. In the first 0.86 My, the dynamical evolution 
of the cloud C displays similar features to that of cloud B. A type A rim firstly forms after 
t = 0.365 My, and then a type B rimmed morphology develops when t = 0.86 My, as shown 
in the corresponding panel in Figure [H The relevant panel in Figure [10] shows the formation 
of an isothermal core with a central density slightly lower than 10^ cm~^ and a thin layer of 
PDR surrounding the core. The average temperature is 35 K in the core and Tj = 1000 
in the photo-ionised bright layer. The thermal pressure P/k in the different regions is ~ 10^. 
The mass included within the isothermal core (PDR included) is 11.1 M0, which is higher 
than the 5 Mq estimated by Equation (fT9|) Therefore, the core would collapse finally. 

However, the gas particles at the two sides of the rear hemisphere continue their move- 
ment to z axis to form a longer tail as shown in the bottom-middle panels of Figures [9] and 
[TT] so that the morphology of the globule becomes a type C BRC at t = 1.034 My. 

The last two panels in Figure [10] shows that the gas particles in the core region collapse 
to the G centre, so the central density reaches to 10^ cm~^, which means the beginning of 
the triggered star formation at t = 1.12 My. 

In table 1. the final morphologies of the simulated clouds are list in the 7th column and 
the 8th and 9th columns are the isothermal Core Formation Time (C.F.T.) and the Core 
Collapse Time(C.C.T.) respectively. 

In order to find the mechanism responsible for the differences in the morphological 
evolution of cloud A, B and C, a further investigation on the kinematics of the shocked gas 
in the front hemisphere is carried in the following part of the paper. 

IV: The Kinematics of morphological evolution 

Figure [5] shows that the direction of the total velocity Vr of the shocked gas particles 
within two 'ear-like' structures determines their converging points on the z axis, hence the 
length of the tail or say it is the magnitudes of Vc and Vq which determine the morphological 
evolution of a BRC. 

The ratio of the gravitational accelerations of cloud A, B and C on the gas particles 
in the 'ear-like' structure is : : ac = 1 • 0.14 : 0.07. Hence the relevant velocity 
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components obey the inequality Vg{A) > Vg{B) > Vg{C). On the other hand, as stated 
in Equation(20), we have Vc{9) oc Vc{0) ~ n"^/^ for fixed Jq and 9, hence for the particles 
in the the 'ear-like' structures in cloud A, B and C, the velocity component Vc follows 
Vc{A) < Vc{B) < Vc{C). Now it is not difficult to understand that the angle 5 of a gas 
particle in the 'ear-like' structure in cloud A, B, C obeys the inequality: S{A) > S{B) > S{C). 
The distance d between the converging point C and the new gravitational centre G in clouds 
A, B and C obeys (Ia < < dc, therefore the three clouds evolved to type A, B, and C 
BRCs respectively. 

On the other hand, the escape velocity of a particle at the surface of a cloud is Ves = 



The ratio of the escaping velocities of three clouds is VA{es) '■ VB(es) '■ vc{es) = 1 : 0.37 : 
0.26. The easiness for the mass particles at the surface of the cloud to be photo-evaporated 
is in the order of cloud C, B, and A, which explains why the masses left in the isothermal 
cores are 11.1, 13.5 and 15 Mq for clouds C, B, and A respectively. 

V: Further increment of the radius 

We further increase the initial radius of the molecular cloud so that the cloud is initially 
much less gravitationally bounded than cloud C. The basic features in the evolutionary 
sequences of the molecular clouds are not much different from that of the cloud C if < 4.7 
pc except that both the mass accumulated in the head of the quasi-stable globules and the 
central density decrease with the initial radius increasing from 1.9 to 4.7 pc. 

However when the initial radius increases to R = 4.7 pc (the cloud is called cloud D in 
Table 1), the initial gravitational force to the surrounding gas particles is only about 1.1% 
of that of cloud A. Cloud D quickly evolves from type A to B then into a loose type C 
BRC - a cometary globule with a slightly compressed core of a central density of a few 10^ 
cm~^ when t = 0.8 My. Shortly after a weak shock passed the centre of the head, the whole 
globule undergoes a quick re-expansion and is totally photoevaporated after 1 My in a so 



called 'ionising flush' (ILefloch &: Lazarefj|l994l ). 



The above set of simulations prove that changing the initial size of the cloud hence the 
initial self-gravity could dramatically change the cloud's evolutionary sequence and the flnal 
morphology. The result means the initial self-gravity of the cloud plays a dominant role in 
the morphological evolutions of BRCs. However, changing of the initial self-gravity of the 
cloud by changing its initial radius creates a second effect at the same time - the alteration 
of its initial ionisation/recombination state, as shown in Table 1. Hence our confldence in 
the role played by the self-gravity in determining the morphological evolution of BRCs needs 
further assurance. The second set of simulations therefore is designed to change the initial 
density of the cloud by changing its initial mass so that the change in the initial self-gravity 
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will not bring drastic change in the ionisation/recombination state as shown in Table 2. 
Hence the role of the self-gravity in the morphological evolution of the BRCs can be further 
explored. 

4.i.2. Clouds with a fixed radius of 0.5pc 

Table 2. lists the parameters of the second set clouds called A, A', A" respectively with 
their initial locations in ( A, F ) parameter space and the information on their morphological 
evolutions. 

For cloud A, the evolutionary features has been fully discussed before, and are re-hsted 
in Table 2 for comparison. With the decrement of the initial mass of the cloud, the density 
of the cloud decreases accordingly, the morphological evolution follows a similar pattern to 
that in the first set of clouds. The dynamical evolutions of the formation of type A, B and 
C BRCs exhibit similar sequences to that in the first set clouds. Therefore in order to save 
space, we don't present the sequential figures on their dynamical evolutions as we did for set 
one clouds, but only provide a summary on the relevant key features which are partly listed 
in the Table 2. 

When the initial mass of a cloud decreased to 12 Mq, cloud A' evolved from type A to 
B morphology in a very similar way to cloud B in set one simulation. An isothermal core 
formed in 0.42 My, and collapsed in 0.8 My, with the central density up to lO'^cm"^ and a 
final mass of 6.9 Mq in a collapsed core of 0.1 pc in radius. When the initial mass further 
decreased to 8 Mq, cloud A" evolved through a sequence of from type A to B then to C 
morphologies, with formation of an isothermal core at t — 0.53 My, which collapsed in 1.5 
My with a central density of lO^cm"^ and a mass of 4.5 M© in a core of radius of 0.1 pc. 

Finally when the initial cloud mass further decreased to 1 Mq, cloud A'" evolved quickly 
through type A and B morphologies with a condensed core of a central density of 10^ cm~^. 
After evolving to a cometary globule morphology at t = 0.3 My, it re-expands and then is 
totally evaporated when t = 0.9 My. 

In the second set simulations, the initial recombination parameter F is same for cloud 
A, a' and A" because they have same radius. The ratio of the initial ionisation parameters 
A{A) : A{A') : A{A") = 1 : 2.85 : 4, i.e., cloud A" has the highest value of A which is 4 
times the lowest value for cloud A. On the other hand, the ratio of the initial gravitational 
energies of the three clouds is Qa '■ ^A' '■ ^A" = 1 : 0.12 : 0.052. The highest value of the 
initial self-gravitational energy of cloud A is about 20 times of the lowest value for cloud A". 
It is obvious the difference in the initial self-gravitational energies among the three clouds 
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is much more distinctive than the difference in their initial ionisation states. Therefore, 
we can safely conclude that the different morphological evolutions of the three clouds are 
caused by their different initial self-gravitational states. The role of the initial self-gravity 
in determining the morphological evolutions of the BRCs is further confirmed. 



4.1.3. The Possible modification of the fA, parameter space 

From the morphological evolution paths of the 7 clouds, we see that for clouds starting 
from region V such as cloud A and A', their dynamical evolutions are driven initially by a 
weak-R ionisation front which quickly changes into a D-type ionisation front, just as what 
Lefloch described. 

However for clouds starting from region III such as cloud B,C,and A", they are not 
at all entirely photo-ion i sed by an R-weak front and permeated by an ionisation flash like 



what iLefloch fc Lazarefj (jl994j ) described, their dynamical evolutions are driven by a similar 



mechanism to that of cloud A and A'. 

For clouds A " and D, which also start from region III but with a high initial ionisation 
state A > 23, its dynamical evolution is driven by a weak R-front and t he whole structure 



get pe rmeated by an ionisation flash and evaporated quickly just as what lLefloch fc Lazarefl 



(I1994J ) described for region III cloud. 



Hence it is clear that in term of the dynamical classification of (A, F) parameter space, 
the boundary for region HI has been pushed greatly from A < ^/lO to A < 23 due to the 
inclusion of the self-gravity of the cloud. 



4.1.4. Zero self-gravity and surface instability of the cloud 

In order to further support our argument on the role of self-gravity of the cloud in the 
morphologic evolution process, we take off the self-gravity of the three clouds A, B and C 
and re-run the simulation s. The result reveals very similar evolutionary sequence to that of 



Lefioch fc Lazarefj (119941 ). The evolutionary sequence passes type A to B morphologies and 
then form a type C BRC. Then it re-expands and finally totally evaporated between 1 to 2 
My. 

To some stage of the evolution of the molecular cloud, the ionisation flux onto the front 
surface of the cloud is perturbed somehow due to the recombinations of the electrons onto 
the ionic atoms, which causes instability of the front surface of the cloud as described in 
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the lefloch's simulations (ILefloch &: Lazarefilll994l ). We found from our simulations that the 



scale of the instability is closely related to the initial ionisation state of the cloud. For all 
of the 7 clouds we simulated, except cloud D, their initial ionisation state are not very high 
when compared to that of the cloud D, very small scale perturbations (with wavelength 
A < 0.2R << R)) at the front surface appeared, as seen from the up-middle panels of Fig.2, 
6 and 9, but they are quickly smoothed out. 

For cloud D, due to its extremely high initial ionisation state, the surrounding of the 
front surface quickly becomes recombination dominant region, so large scale perturbations 
(with wavelength 0.2R < X < R) grow at t = 0.06 My from the front surface of the cloud 
which can be seen from the up- left panel of Figure [T2j However the large scale perturbation 
is stable because it is gradually smoothed out with the core further being compressed by the 
isothermal shock. At t = 0.8 My, the front surface becomes very smoothed. We think it is 
the gradually increased gravity from the core underneath the front surface which suppressed 
further growth of the surface instability. 

We conclude that with inclusion of the self-gravity in the model, both small and large 
scale surface instabilities caused by the perturbation to the ionisation flux are stable and will 
be smoothed out over the further evolution of the c loud. Our conclusion is in agreement with 



what revealed by iKessel-Deynet fc BurkertI (120031 ) in Case B in their paper on the stability 



of the front surface of the cloud. 

We have shown that the self-gravity indeed plays an important role in the morphological 
evolution of the molecular clouds. For any change in the initial conditions of a cloud, if 
it greatly results in a change in its initial self-gravitational state, the cloud will alter its 
morphological evolution. 



4.2. Initially Thermally Supported Clouds 

In order to see how the three clouds would evolve if they are only initially thermally 
supported, i.e., if the cloud A, B and C start their evolutions without the initial turbulence. 
We re-run the simulations for cloud A, B and C respectively by setting a zero initial velocity 
field to each of them. Since the general features of the morphological evolutions of the 
three molecular clouds are qualitatively similar to that presented in last section, we mainly 
describe the different features in their evolutions. 

Cloud A evolved into a type A BRC with a condensed core of a central density of a few 
10®cm~'^ within 0.37 My. Cloud B evolved into a type B BRC with a condensed core of a 
central density lO^cm"^ within 0.74 My. Cloud C reached a type C BRC structure with a 
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core density of lO^cm"^ within 0.98 My. Compared with the corresponding results in Section 
4.1.1, the clouds without initial turbulence take less time to collapse. Further comparing the 
velocity field evolutions of the two group molecular clouds, we find the initially turbulent 
field takes 0.02 My in cloud A, 0.07 My in cloud B and 0.14 My in cloud C to dissipate, 
which shows the initially turbulent clouds is just delayed to form a highly condensed core. 
Therefore we conclude that the initial turbulence does not effect a cloud's final rnqrphol ogy 
type. The results is consistent with that of an investigation by iNelson fc Langen (119991 ) on 
the dynamical evolutions of Bok globules affected by FUV radiation. 



4.3. Effect of changing initial density profiles 

All of the above simulations are based on an initially gravitationally stable and uniform 
cloud. Although large degrees of central condensation should not be expected in the ini- 
tially gravitationally stable molecular cloud, we are still interested in examining the effect 
of the initial central condensation of a molecular cloud on its morphological evolution. The 
following initial density distribution is used, 

_ f '^o (0<r <ro) 
^ ^ I rio(^)^ (ro < r < R) 

which represents for a condensed uniform core with a radius of ro and density of no sur- 
rounded by an envelop of a density profile proportional to r~^. For a given ro, R and M, no 
can be found from the equation AwmH n{r)r'^dr = M. 

When we start with ro = 0.5-R, the evolutionary sequences of clouds A, B and C with 
an initial core-halo mass distribution do not show significant change from that of initially 
uniform clouds. The basic features of type A, B and C are kept, except that the main body of 
the formed BRC structure are all smaller than that formed in initially uniform clouds. Then 
we decrease the radius of the central condensation to ro = 0.2R, the morphological evolution 
show obvious changes. Cloud A ends its evolution in type B as shown in Figure [13], and cloud 
B in type C as shown in Figure [lH while cloud C developed into a very narrow shaped type 
C at a much early time ( t = 0.67 My ) as shown in Figure [T5l Their morphologies are all 
elongated than that of their corresponding uniform clouds. This morphological change can 
be well explained by the illustration in Figure ([5]), the gas particles in the halo has a lower 
density therefore a higher Vc component, which makes them converge to a farther point C 
from point G on the z axis. 

The mass left at the collapsed core of 0.15 pc in cloud A, B, and C are 10.5, 7.2, 5.6 Mq 
respectively, which are less than that in their corresponding initially uniform clouds. This is 
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because the density in the initial halo of a cloud is lower than that in the same part of the 
initially uniform cloud hence the halo in the front hemisphere is at a higher initial ionisation 
state, in which particles are easier to be ionised and then to be photo-evaporated. 

In summary, the initial mass condensation state does affect a cloud's morphological 
evolutions because it changes the initial self-gravitational state of a cloud. 

4.4. Effects of the boundary conditions 

As stated in Section 2.5, we apply a constant pressure boundary condition, which re- 
sembles the physical conditions of a very hot HII region. In order to examine the effect 
of the boundary condition on the dynamical evolution of a cloud we decrease the value of 
the surrounding pressure Pext/k to mimic a less violent environment. For each of three 
clouds in the first set, we re-run the simulation for two different values of the external pres- 
sure, i) Pp^tlk = 10^ Kc m~^ which corresponds to a warm interstellar medium environment 



(INelson fc Langerlll999l ) and ii) Pext/k = which describes a vacuum environment so that 



the ionised gas is able to freely expand from the surface of the cloud. 

It is found that with decreasing external pressure, although cloud A still evolves to a type 
A BRCs under the boundary condition i) and ii), the curvature of the front surface increases 
and the shape of the final condensed core becomes elongated. Cloud B with boundary 
condition i) evolves to a slightly elongated type B BRC but evolves to a type C BRC at 
t = 1.6 My under the boundary condition ii), as shown in Figure [161 Cloud C forms a type 
C BRC in both cases, but the length of the tail increases with the decrease of the external 
pressure. It is obvious that the lower pressure environment favours formation of morphology 
with higher curvature, or high pressure environment favours formation of type A BRCs. The 



concl usion is consistent with that of zero velocity boundary condition (ILefloch fc Lazarefl 



19941 ). 



4.5. The effect of initial temperature of molecular clouds 



According to ob servations, the temperatures of molecular clouds are generally in the 
range of 10 — 100 K (jBertoldilll989l ). We then carried out two more sets of simulations for 
each of the three molecular clouds A, B and C with the initial temperatures Ti = 20 K and 
T2 = 100 K respectively, while all of the others physical conditions are kept the same as that 
in Section 4.1.1. 



The evolutionary sequence of each simulation is qualitatively similar to its corresponding 
cloud with an initial temperature of 60 K. 
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When the initial temperature Ti = 20 K, the evolutionary sequences of cloud A, B and 
C basically exhibit no significant differences from that described in Section 4.1.1. Cloud A 
evolves to a quasi-stable type A rimmed morphology with a condensed core at the head of 
its front hemisphere. The central density of the core reaches a value of 10^ cm~^ in 0.37 My 
when it finally collapses. Cloud B evolves from type A to B morphology with a condensed 
core at the head of its front hemisphere in 0.6 My and the core reaches a central density of 
10^° cm~'^ at t = 0.79 My. Similarly cloud C evolves to a type C BRC after passing type 
A and B morphologies. It collapse at t = 1.01 My when the central density of the core is 
2 X 10^ cm-^ 

When the initial temperature rises to 100 K, cloud A still evolves to a type A rimmed 
morphology having a condensed core which collapses after the central density is up to a few 
10^ cm~^ in the front hemisphere at t = 0.42 My. Cloud B evolves to a type B rimmed 
morphology and the condensed core collapse with the central density being up to 10^ cm~'^ 
at t = 0.93 My. Cloud C evolves to a type C rimmed BRC plus a small core of centre density 
of 10^ cm-3 in 1.32 My. 

In summary, the variation of the initial temperature of molecular clouds between 10 - 100 
K does not significantly change the morphological evolution of a cloud, which is physically 
understandable sin ce the thermal pressure of a cloud is not important in the dynamical 



evolution in BRCs flBertoldilll989l : lLefloch fc LazarefBll994l : IWhite et al.lll999f ). Higher initial 



temperature in a cloud delays its collapsing time, for the initial thermal motions of gas 
particles dissipate over a period of time. The morphological evolution of a cloud is mainly 
decided by the dynamical behaviour of the gas particles within the two sides of the cloud 
{6 ~ 7r/2), which is largely dependent on the strength of the ionising fiux and the self-gravity 
of the cloud as shown in Section 4.1. 

However, the initial temperature of a cloud does affect the amount of mass condensed 
in the core of a BRC, for which we will have a detailed discussion when we study the UV 
radiation triggered star formation in BRCs in our next paper. 



4.6. The effect of the strength of ionising flux 

We also carried out a set of simulations on the influence of the strength of an ionising 
flux on the evolutionary process of BRCs. The change in the strength of ionising flux incident 
on the front surface of a molecular cloud corresponds to two different cases: i) a cloud can be 
at different distances from an ionising star; ii) an ionising star in different types radiates with 
different flux strength. This set of simulations is to see whether the molecular clouds within 
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one giant cloud cluster would exhibit essential difference in their morphological evolutions 
due to their different distances from the central star. 

Since the ionising flux J(Lyman) oc 1/-D^, where D is the distance of a cloud from the 
ionising star. For each of cloud A, B and C, we re-run the simulations twice by setting two 
different J(Lyman) flux values. If we assume cloud A, B and C in Section 4.1.1 are in a 
distance of Dq from the ionising star, then Ji (Lyman) = Jq (Lyman) /4 means we move the 
star twice distance away, i.e., D = 2Dq. If we set J2 (Lyman) = 4 x Jo (Lyman) which means 
the original cloud is moved to half distance closer, i.e, D = Do/2. The spatial range of 0.5 - 
2 Dq from the ionising star should be large enough to include most molecular clouds in one 
cloud cluster. 

The simulation results revealed that for all of the three molecular cloud A, B, and C, 
their morphological evolutions are not very sensitive to the strength of the ionising flux in 
the above specified range and they follow very similar evolutionary sequences to that when 
ionising flux is Jo{Lyman) . 

The result looks a bit puzzling at the first glance, since the velocity component Vc 
of the ionised gas in the top layer of the cloud increases with the ionising flux, as shown 
in Figure O so that we should expect the angle 6 decreases with the ionising flux, which 
will favour the formation of type C rim morphology. However, if we consider the effect 
of the enhanced shock effect due to increased ionising flux, the compression to the core in 
the front hemisphere will be stronger. Therefore the core will exert a stronger attractive 
force on its peripheral gas particles which then results in a higher Vg so that angle 6 of 
the velocity remains basically unchanged although the magnitude of the total velocity may 
increase, which would decrease the time needed for the gas particles in the ear-like structures 
to collapse onto the symmetrical axis. The same consideration can explain the insensitivity 
of the morphological evolution of the cloud to the decreased radiation flux and the time 
needed for the gas particles to collapse onto the symmetrical axis is found to be elongated. 

4.7. Confrontation with observations 

The simulations on a set of molecular clouds reveal that the morphological evolution of 
a molecular cloud under the effect of the UV radiation does not necessarily go through all of 
the three morphologies and then end up at cometary type C morphology as described by the 
previous models. Depending on its initial gravitational state, a cloud could evolve to any one 
of the three type BRCs. The diversity of the evolutionary sequence of BRCs revealed from 
our modelling provides a good prospect for us to explore the possible explanations for the 
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questions raised from BRCs' observations and presented in the introduction of this paper. 

The modelhng results for the molecular cloud A, B, C under the same initial temper- 
ature and ionising flux tell us that the evolutionary paths of all molecular clouds would 
firstly pass the type A morphology. However some of molecular clouds would only evolve 
to a quasi-equilibrium type A morphology plus a highly condensed core at the head of its 
front hemisphere, if their initial self-gravitational energy is high enough. Therefore we shall 
statistically be able to observe many more type A BRCs than type B or C BRCs in a random 
observation. 

Our modelling results also shows that the UV radiation triggered collapse of a cloud onto 
the centre core in its front hemisphere could occur with any one of the three morphologies, 
so it is not surprising at all that the UV radiation triggered star formation occurring in type 
A BRCs could be observed. 

Finally, we come to understand the spacial sequen ce of type A-B-C BRCs with their 



distance from the centre star in the Ori OB 1 association (jOgura fc Sugitanilll998[ ). We might 
think this morphological distribution of BRCs as the result of decreasing illuminating flux 
with distance. However as shown in the argument in Section 4.6, the morphological evolution 
of the cloud is insensitive to the change of the UV radiation flux in the range discussed. Based 
on the results derived from the simulations, we can understand this observation from two 
aspects. 

First, we know that where a new star is formed inside a giant molecular cloud as a 
result of local molecular gas collapse onto a central point, the density of the surrounding 
gas can be reasonably assumed to decrease with the distance from the cent ral star. It is 



also known that giant molecular clouds are very clump (ITatematsu et allll99ll ). Considering 
the above two points we can assume that the densities of molecular clumps surrounding a 
new star decrease with their distance D from the central star. If we further assume that 
the clumps in a giant molecular cloud are of similar radius i?, then the initial gravitational 
potential energy of a clump decrease with their distance from the centre of the star for the 
gravitational potential energy can be written as ~ M^/i? ~ R^n"^, where n is the initial 
number density of the clumps. 

The molecular clumps closest to the central star could possibly evolve to a quasi- 
equilibrium type A BRC because they have highest initial self-gravitational potential en- 
ergies, while for those clumps farther away, they evolve from type A to a quasi-equilibrium 
type B BRC because of their moderate initial self-gravitational potential energies. For those 
farthest away from the central star, they evolve to pass type A and B morphologies, then to a 
cometary type C morphology due to their lowest initial self-gravitational potential energies. 
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Secondly, we already showed in Section 4.4 that at the borders of HII region, type A 
BRCs are favoured. 

Therefore the spatial sequence of type A - B - C BRCs with their distance to the central 
star can be seen as a manifestation of the mass density distribution of the very clumpy giant 
molecular clouds. 

5. Conclusions 

A three dimensional and comprehensive RDI model has been developed with the SPH 
technique, which self-consistently includes the hydrodynamical evolution, the self-gravity of 
the cloud, the energy evolution, radiation transferring and a basic chemical network. T 
The investigation on the morphological evolutions of BRCs are carried out by employing 
the newly developed model in order to understand a series of puzzles coming from BRCs 
observations. It is found that the morphological evolution of a molecular cloud is much more 
sensitive to its initial self-gravitational state than to other physical conditions such as the 
initial temperature between 20 and 100 K, or the ionising radiation flux within one order 
of magnitude higher or lower than Jo{Ljma.n). Our modelling results have revealed that 
there are three different evolutionary prospects when a molecular cloud is under the effects 
of ionisation radiations. 

For a molecular cloud of an initial mass 35 Mq and initial temperature 60 K under 
the effects of UV and FUV radiation fields specified above, the simulation results revealed: 
a) if its initial radius is 0.5 pc, it will evolve to a quasi-stable type A BRC; b) when its 
initial radius increases to 1.3 pc, its evolution will firstly pass a type A morphology and 
continue its journey to a quasi-stable type B rimmed morphology; c) when the initial radius 
further increases to 1.9 pc, its evolutionary process follows a sequence of type A— >-B— >-C 
morphologies, d) Further increasing the cloud's radius to 4.7 pc, the evolution of the cloud 
also follows the sequence of type A — >B^C, but the finally formed type C BRC will be 
totally photoevaporated. 

Further simulation results also show that the initial central condensation of the mass in a 
cloud could change its morphological evolution in term of its initially uniform correspondence. 
At the border of HII regions, the formation of type A BRCs are favoured. When a molecular 
cloud is set to a zero initial self-gravitational potential energy or with a low initial density, 
its morphological evolution sequence is consistent with that by Lefioch's model which is a 
2-D model and didn't include the self-gravity of the cloud. 

With the knowledge obtained from the modelling, we are now able to confront with the 
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questions raised by observations. Since all of the clouds evolve to pass a type A morphology 
and some of the clouds with enough high initial gravitational potential energy would just 
evolve to a quasi-stable type A morphology, statistically many more type A BRCs should be 
observed over a random observation. For those BRCs who evolve to a quasi-stable type A 
rimmed morphology, if the core mass contained in their front hemisphere is high enough for 
it to collapse, new star could be triggered to form. Star forming signal should be observed. 
Furthermore if the initial densities of the molecular clumps in Ori OB 1 associations decrease 
with their distances from the centre of the ionising star, a spacial distribution of type A-B-C 
BRCs with distance from the central star could exist. 

It is obvious that the inclusion of the self-gravitation in our model revealed a diversity 
of the evolutionary sequence of a molecular cloud under the effects of UV radiations. The 
satisfactory explanations based on the modelling results to the observational puzzles state 
that the self-gravity of a molecular cloud indeed plays an important role in controlling 
the evolutionary destiny of a molecular cloud and should not be neglected especially when 
we try to achieve a whole picture of the dynamical evolutions of molecular clouds. This is 
because the radiation induced shock dramatically enhances the effect of the self-gravity of the 
clouds. Consequently, some modifications on the region classification in Lefioch's parameter 
space have been confidently derived. The boundary between the D-type ionisation front and 
'ionising flush' is extended from 5 — y/TO to 23 as a result of the inclusion of the self-gravity 
of the molecular cloud. 
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Fig. 1. — The initial geometry of a molecular cloud and the configuration of the radiation 
fields. The heavy arrow lines above the upper hemisphere represent the sum of the strong 
UV (Lyman continuum) fiux from the massive stars above the molecular cloud and the 
isotropic FUV fiux from the background star radiation and the light arrow lines below the 
lower hemisphere represent the isotropic FUV radiation {hu < 13.6 eV) fiux only. 
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Fig. 2. — The evolutionary snapshots of the number density of the simulated molecular cloud 
A at the cross section y = from t = 700 years to t = 0.394 My. The top grey bar shows the 
density scale in logarithm for the upper row snapshots and the bottom bar for the second 
row snapshots. 
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Fig. 3. — The distributions of the number density (sohd hne) and ionisation ratio x (dot- 
dashe hne) of the simulated molecular cloud A along the symmetrical axis z from t — 700 
years to t = 0.394 My. 
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Fig. 4. — The evolution of the velocity field for the same cloud as in Figure [2] over a period 
of 0.394 My. The length of the arrow indicates the magnitude of velocity and the maximum 
magnitude of the velocities at different stages of the evolution have been written inside the 
corresponding panels. 
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Fig. 5. — An illustration about the effect of the gravity of the main body of the cloud on the 
velocity of the shocked gas particle in the 'ear-like' structures at the two sides of the cloud. 
The central solid lines are their structure symmetrical line z. 
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Fig. 6. — The evolutionary snapshots of the number density of the simulated molecular cloud 
B at the cross section y = from t = 1800 years to t = 0.81 My. The top grey bar shows the 
density scale in logarithm for the upper row snapshots and the bottom bar for the second 
row snapshots. 
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Fig. 7. — The distributions of the number density (sohd hne) and ionisation ratio x (dot- 
dashe Hne) of the simulated molecular cloud B along the symmetrical axis z from t — 1800 
years to t = 0.81 My 
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Fig. 8. — The evolution of the velocity field for the same cloud as in Figure [6] over a period 
of time of 0.81 My. The length of the arrow indicates the magnitude of velocity and the 
maximum magnitude of the velocities at different stages of the evolution have been written 
inside the corresponding panels. 
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Fig. 9. — The evolutionary snapshots of the number density of the simulated molecular cloud 
C at the cross section y = over a period of 1.12 My. The top grey bar shows the density 
scale in logarithm for the upper row snapshots and the bottom bar for the snapshots in the 
bottom row. 
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Fig. 10. — The distributions of the number density (sohd hne) and ionisation ratio x (dot- 
dashe Une) of the simulated molecular cloud C along the symmetrical axis x — over a 
period of 1.12 My 
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Fig. 11. — The evolution of the velocity field for the same cloud as in Figure [9] over the 
same period of time. The length of the arrows indicates the magnitude of velocity and the 
longest arrow in each panel represents the maximum velocity in the corresponding velocity 
field distribution at different stages of the evolution. 
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Fig. 12. — The growth and smooth processes of the large-scale surface instability over the 
evolutionary process in cloud D. 
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Fig. 13. — The evolutionary snapshots of the number density of the simulated molecular 
cloud A with an initial core-halo distribution at the cross section y = 0. The top grey bar 
shows the density scale in logarithm for the upper row snapshots and the bottom bar for the 
second row snapshots. 
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Fig. 14. — The evolutionary snapshots of the number density of the simulated molecular 
cloud B with an initial core-halo distributionat the cross section y = 0. The top grey bar 
shows the density scale in logarithm for the upper row snapshots and the bottom bar for the 
second row snapshots. 
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Fig. 15. — The evolutionary snapshots of the number density of the simulated molecular 
cloud C with an initial core-halo distributionat the cross section y = 0. The top grey bar 
shows the density scale in logarithm for the upper row snapshots and the bottom bar for the 
second row snapshots. 
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Fig. 16. — The evolutionary snapshots of the number density of the simulated molecular 
cloud B with a zero external pressure boundary at the cross section y — 0. The top grey bar 
shows the density scale in logarithm for the upper row snapshots and the bottom bar for the 
snapshots in the bottom row. 
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^Morphology type 

^Core formation time in My 

^Core collapse time in My 



Table 1. Parameters of set one clouds 



Cloud 


A 


r 


no (cm ^) 


R(pc) 


Region 


M type^ 


CFT2 


CCT3 


A 


0.68 


55 


2672 


0.5 


V 


A 


0.35 


0.39 


B 


7.3 


89.1 


152 


1.3 


III 


B 


0.63 


0.6 


C 


18 


108 


49 


1.9 


III 


C 


0.86 


1.12 


D 


182 


169 


3.2 


4.7 


III 


c 


1.0 


N/A 
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They have the same meanings as that in Table 1. 



Table 2. Parameters of second set clouds 



Cloud 


A 


r 


no (cm ^) 


M(Mq ) 


Region 


M type^ 


CFT2 


CCT3 


A 


0.68 


55 


2672 


35 


V 


A 


0.35 


0.39 


A' 


1.94 


55 


916 


12 


V 


B 


0.45 


0.55 


A" 


2.9 


55 


610 


8 


III 


C 


0.48 


0.72 


A'" 


23 


55 


76 


1 


III 


C 


0.25 


N/A 



